Determining spatial structures of ion crystals by simulated annealing method
Wu Wen-Bo1, 2, Wu Chun-Wang1, 2, Li Jian1, 2, Ou Bao-Quan1, 2, Xie Yi1, 2, Wu Wei1, 2, †, Chen Ping-Xing1, 2, ‡
College of Science, National University of Defense Technology, Changsha 410073 , China
Interdisciplinary Center for Quantum Information, National University of Defense Technology, Changsha 410073 , China

 

† Corresponding author. E-mail: weiwu@nudt.edu.cn pxchen@nudt.edu.cn

Project supported by the National Basic Research Program of China (Grant No. 2016YFA0301903), the National Natural Science Foundation of China (Grant Nos. 11304387, 11174370, 61632021, 61205108, and 11305262), and the Research Plan Project of National University of Defense Technology (Grant No. ZK16-03-04).

Abstract

Calculating the spatial structures of ion crystals is important in ion-trapped quantum computation. Here we demonstrate that the simulated annealing method is a powerful tool to evaluate the structures of ion crystals. By calculating equilibrium positions of 10 ions under harmonic potential and those of 120 ions under anharmonic potential, both with the standard procedure and simulated annealing method, we find that the standard procedure to evaluate spatial structures is complicated and may be inefficient in some cases, and that the simulated annealing method is more favorable.

1. Introduction

The trapped ions system serves as one of the most promising systems for realizing an operational quantum computer.[1] Among all kinds of configurations of the ion traps, the surface-electrode ion trap is distinctive due to its versatility, such as the high scalability and the good stability of the trap.[2,3] Recently, the large-scale ion traps which have been developed from the quadrupole trap attract great attention.[4,5] These traps have quantum charge-coupled device (QCCD) architectures, in which ions are cooled in one region but stored and used in other regions. In a particular region, using a high frequency harmonic potential along the axial direction allows for the ion string experiencing a second-order phase transition from a linear to a zig–zag structure.[68] In fact, after a structural phase transition, the trapped ions will always be in certain equilibrium positions with the minimum total potential energy.[911] Between two particular regions, ions can be transported adiabatically, it is a slow transport process in which ions are in equilibrium positions at every moment. So when we discuss the trapping stability, structural phase transition and transportation of the ions, we need to calculate the equilibrium positions of the ion crystals.

There is a standard procedure of calculating the equilibrium positions of the ion crystals, its details are given by James.[12] Many groups have adopted the standard procedure no matter the axial quasistatic electric field is harmonic[12,13] or anharmonic.[14,15] However, the standard procedure may be hard to use under some particular axial potentials, and it may not find the right results. Thus to find a better algorithm is necessary. Simulated annealing (SA) method is a heuristic algorithm, it simulates the heated metal’s cooling process which is somewhat analogous to cooling ions, and it is powerful to solve this kind of optimization problems.

In this article, the SA method is introduced for computing the equilibrium positions of the ions. Then we compare the standard procedure and the SA method by calculating equilibrium positions of 10 ions under harmonic potential, we also compare the two methods when calculating equilibrium positions of 120 ions under anharmonic potential. We find the SA method is more powerful than the standard procedure to evaluate the spatial structures of ion crystals.

2. Calculating equilibrium positions of ion crystals

A linear ion string structure in Paul trap is the simplest type of ion Coulomb crystals.[16] It arises mainly in the case of the Coulomb interaction between every two ions, and the external electric field, which is generated from a radial radio frequency (RF) and a much weaker axial electrostatic field along the trap. We obtain the structure by calculating the equilibrium positions of the ions under the condition that the system’s potential energy is minimized. We consider a string of N ions each with charge q. The position of the m-th ion is denoted by , and the ions are numbered from left to right. Hence the potential energy of the system is given by where V denotes the electrostatic potential along the trap axis, and is the permittivity of free space. By means of researching the global minimum of the potential energy U, the equilibrium positions of the ions can be found.

2.1. Standard procedure and simulated annealing method
2.1.1. Standard procedure

Briefly speaking, the standard procedure of minimizing the potential energy has three steps: (i) taking the partial derivative of potential energy for each position variable, (ii) making them equal to zero, and (iii) sovling these partial differential equations. Whether the axial electrostatic potential V is in the form of harmonic potential[12,13] or anharmonic potential,[14,15] it is summarized as a process of solving the following set of partial differential equations With V denoting a harmonic potential form, these equations can be solved analytically only for and ; for larger N they can only be solved numerically.[12] For , the results have been tabulated by James,[12] and have been confirmed by many experiments.[13] If V is a complex potential, these equations may only be solved numerically. The standard procedure is straightforward to implement, however, it has two problems: i) calculating partial derivatives of the potential energy U can be very complicated, especially when V has a complex, segmented or discrete form; ii) the positions of ions may be found at the local minimum instead of the global minimum of the potential energy U.

2.1.2. Simulated annealing method

To avoid the problems of the standard procedure mentioned above, we introduce the SA method[17,18] to minimize the potential energy. The SA method is an algorithm motivated by the thermodynamical process, which is cooling molten metal to its state with minimum energy slowly.[19,20] This method has the ability to avoid poor local minimum.[21] Besides, it needs not taking the derivative of U.

The SA algorithm is briefly shown in Fig. 1. It can be summarized as follows: (i) Define the formula of potential energy U and the positions of the ion string , set the halting temperature , the halting potential energy , higher bound of positions , lower bound of positions , the maximum step number , and the cooling rate . Then initialize the temperature and the starting positions , and calculate the potential energy U with these parameters. (ii) Make randomly move to its neighbor , which means the position of each ion ( ) will randomly change to a new value ( ), respectively. Then calculate the new potential energy . (iii) If is lower than halting potential energy , we believe is low enough and the new positions represent the ground state of the ion string. (iv) Otherwise, we accept with a probability of , and the step number K increases by one, where . Here the Boltzman constant is included in T. (v) If step number K is smaller than , we turn to the step (ii). (vi) After searching for times, we drop the temperature by . (vii) If T is still higher than , turn to the step (ii) and repeat the circuit.

Fig. 1. The simulated annealing algorithm.

Here we specially discuss the process of escaping from local minimum to global minimum by the SA method. If new potential energy is lower than or equal to U, probability P is greater than or equal to 1, then the new positions are completely accepted; if is higher than U, though are not better than , we still accept with the probability P which is less than 1 but more than 0, and start next searching cycle from , that is the reason why the SA method has a chance to escape from local minimum and finally find the global minimum. As the temperature T decreasing, the probability P is declining, then new positions are more and more difficult to accept, that will make positions converge.

2.2.2. Calculating a linear ion string structure under harmonic potential

In order to show the validity of our SA algorithm, we calculate the equilibrium positions of the ions under a harmonic potential by both the standard procedure and the SA method, and we compare our results with the tabulated results which are from James.[12] The formula of the harmonic potential can be expressed by where M is the mass of an ion, is the axial trapping frequency of the harmonic potential.

We use strings of ions trapped at the frequency . The mass of each ion is kg. In both methods we set the initial positions randomly. Different initial positions only affect the searching time. For the standard procedure, we use the Newton–Raphson method to calculate the partial differential equations (2) with the absolute tolerance of . In the SA method, by experience we set the initial temperature J, halting temperature J (for a high halting temperature, the iterative results may still change at the end of the iteration, however, for this temperature the iterative results kept invariant for quite a while before being outputted), cooling rate , halting potential energy J and the maximum step number . We have calculated and compared the results for . For the sake of brevity, we only display the results for N = 10 here. To compare the results from these two methods with the results listed by James,[12] here we also define the length scale and calculate the dimensionless equilibrium position . The comparison of the results are given in Table 1.

Table 1.

Scaled equilibrium positions of 10 ions under harmonic potential.

.

In the table, (I), (II), and (III) represent the results from James,[12] the standard procedure and the SA method respectively. By comparisons, we find that our SA algorithm is valid.

2.3. Calculating a linear ion string structure under anharmonic potential

A long linear ion string with equidistance benefits large-scale quantum processor,[14] because it has the closely spaced transverse phonon modes and good addressability for lasers. Under the same radial confinement, a linear ion string in an anharmonic trap can contain more ions than that in a harmonic trap before undergoing a structural phase transition. Anharmonic traps allow to store linear strings of equal distance.[22] To show that SA method is also applicable to the problem in anharmonic traps, we calculate the equilibrium positions of an ion string under the anharmonic potential, whose analytic expressions have been found.[22] In the Ref. [22], the author showed that N ions can form a string with equal separation d between every two adjacent ions under the anharmonic potential as follows: where is an offset potential, is the gamma function (see Appendix A), , the dimensionless and .

Now we calculate the equilibrium positions and the separations between every two adjacent ions under this potential. Using the standard procedure, it is hard to find a proper numerical method to solve Eq. (3) because of the complex potential form. We have tried the fixed-point iteration method, but various iteration equations we have constructed do not work for this potential form. When we turn to the Newton–Raphson method or the quasi-Newton method, it is inevitable to evaluate a mazy Jacobian, apart of that, new iterative positions are usually out of receivable range which is limited by the gamma functions in this potential form. Owing to the trouble mentioned above, we can not keep the stability of the standard procedure. In principle, the SA algorithm totally does not depend on the specific form of the potential energy. We assume the neighboring distance m, reset J and , initialize positions randomly, and calculate the equilibrium positions with the SA method. Here we calculate the structure of 120 ions, namely we minimize the potential energy function with 120 parameters. Define , denotes the average of , we obtain the distribution of ions as shown in Fig. 2.

Fig. 2. (color online) The distribution of the ion spacings by calculating 120 positions of the ions and 119 separations of every two adjacent ions. The average of the separations m, which is very close to its assuming value, and the separations’ standard deviation .

As what we have anticipated, the SA method is powerful to find the right equilibrium positions.

3. Discussion and conclusion

As is shown above, the SA method can find the right results, instead of falling into some unknown local minimum, and can also compute the structure without the limitation of potential form by avoiding evaluating the partial derivatives of the potential energies, especially when the potential has a complex, segment or discrete form, which is difficult or incapable for the standard procedure. Expending this method to two-dimensional and three-dimensional potentials, we only need to double and triple the number of variables for one-dimensional potential, respectively. Moreover, the SA method can also optimize other physical problems effectively. Using this method one can simulate the zig–zag structure and the transportation of ions on a gold-on-silica segmented surface-electrode ion trap.[23,24] For two-dimensional and three-dimensional occasions, one can study the phase transition of the ions, like Kibble–Zurek mechanism with the SA method.[25,26] In addition, by precisely determining the structures of ion crystals with the SA method, one can also study the ion micromotion, it can further promote designing of ion traps.

In summary, we have introduced the SA method to compute the equilibrium positions and the structures of linear ion string under both the harmonic and the anharmonic potentials respectively. By comparing the SA method with the standard procedure, we have found that the SA method is more versatile and more reliable for finding equilibrium positions of ions under various potential forms.

Reference
[1] Ladd T D Jelezko F Laflamme R Nakamura Y Monroe C O’Brien J L 2010 Nature 464 08812
[2] Amini J M Uys H Wesenberg J H Seidelin S Britton J Bollinger J J Leibfried D Ospelkaus C VanDevender A P Wineland D J 2010 New J. Phys. 12 033031
[3] Shu G Vittorini G Buikema A Nichols C S Volin C Stick D Brown K R 2014 Phys. Rev. A 89 062308
[4] Kielpinski D Monroe C Wineland D J 2002 Nature 417 00784
[5] Monroe C Kim J 2013 Science 339 1231298
[6] Blümel R Chen J M Peik E Quint W Schleich W Shen Y R Walther H 1988 Nature 334 334309a0
[7] Ejtemaee S 2015 Dynamics of Trapped Ions Near the Linear-Zigzag Structural Phase Transition Ph. D. Dissertation Burnaby Simon Fraser University
[8] Klumpp A Liebchen B Schmelcher P 2016 Phys. Rev. A 380 2644
[9] Diedrich F Peik E Chen J M Quint W Walther H 1987 Phys. Rev. Lett. 59 2931
[10] Wineland D J Bergquist J C Itano W M Bollinger J J Manney C H 1987 Phys. Rev. Lett. 59 2935
[11] Drewsen M Brodersen C Hornekær L Hangst J S Schifffer J P 1998 Phys. Rev. Lett. 81 2878
[12] James D F V 1998 Appl. Phys. B 66 181
[13] Morigi G Fishman S 2004 Phys. Rev. 70 066141
[14] Lin G D Zhu S L Islam R Kim K Chang M S Korenblit S Monroe C Duan L M 2009 Europhys. Lett. 86 60004
[15] Home J P Hanneke D Jost J D Leibfried D Wineland D J 2011 New J. Phys. 13 073026
[16] Thompson R C 2015 Contemp. Phys. 56 63
[17] Kirkpatrick S Gelatt C D Vecchi M P 1983 Science 220 671
[18] Press W H Teukolsky S A Vetterling W T Flannery B P Numerical Recipes 3 Cambridge Cambridge University Press 549 555
[19] Goffe W L Ferrier G D Rogers J 1994 J. Econometrics 60 65
[20] Laarhoven P J M Aarts E H L Lenstra J K 1992 Oper. Res. 40 113
[21] Johnson D S Aragon C R McGeoch L A Schevon C 1989 Oper. Res. 37 865
[22] Johanning M 2016 Appl. Phys. B 122 71
[23] Ou B Q Zhang J Zhang X F Xie Y Chen T Wu C W Wu W Chen P X 2016 Sci. China-Phys. Mech. Astron. 59 123011
[24] Xie Y Zhang X F Ou B Q Chen T Zhang J Wu C W Wu W Chen P X 2017 Phys. Rev. A 95 032341
[25] Pyka K Keller J Partner H L Nigmatullin R Burgermeister T Meier D M Kuhlmann K Retzker A Plenio M B Zurek W H Campo A Mehlstaubler T E 2013 Nat. Commun. 4 2291
[26] Ulm S Roßnagel J Jacob G Degünther C Dawkins S T Poschinger U G Nigmatullin R Retzker A Plenio M B Schmidt-Kaler F Singer K 2013 Nat. Commun. 4 2290